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We analytically calculate the energy, magnetization curves (B(H)), and elasticity of skyrmion 
flux lattices in p-wave superconductors near the lower critical field H c i, and use these results with 
the Lindemann criterion to predict their melting curve. In striking contrast to vortex flux lattices, 
which always melt at an external field H > Hci, skyrmion flux lattices never melt near H c \. This 
provides a simple and unambiguous test for the presence of skyrmions. 
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The topological excitations known as skyrmions have 
been proposed to have many applications. After they 
were introduced in a nuclear physics context by Skyrme 
variations of this concept have been shown or pro- 
posed to be important in superfluid 3 He 0, in quantum 
Hall systems 0], in p-wave superconductors Q, and in 
metallic magnets [Bj]. Given this widespread predicted 
occurrence, it is desirable to find a simple experimental 
signature that can serve as a smoking gun indicating the 
presence of skyrmions. In this Letter we show that, in 
p-wave superconductors subject to an external magnetic 
field, the structure of the phase diagram provides an un- 
ambiguous test for the presence of skyrmions Q . 

In superconductors, skyrmions compete with another, 
and more well-known, species of topological excitations, 
viz., vortices. In type-II superconductors, the latter are 
induced by an external magnetic field H and form the 
famous Abrikosov flux lattice in a field range H c i < H < 
H c2 . It is known both theoretically 0, @, El EU and ex- 
perimentally [IH that these flux lattices can melt. 
They do so near both H c i and H C 2, where the elastic con- 
stants of the flux lattice vanish, which in clean systems 
makes the root-mean-square positional thermal fluctua- 
tions \J (|u(a;)| 2 ) grow without bound as these fields are 
approached. 

Vortices involve a singular field configuration at their 
cores and can occur for any symmetry of the supercon- 
ducting order parameter. In p-wave superconductors the 
vector character of the superconducting order parame- 
ter also allows for skyrmions, which, in contrast to vor- 
tices, are non-singular topological defects. Like vortices, 
each skyrmion carries a quantized magnetic flux, and 
for strongly type-II p-wave superconductors in a field 
H c i < H < H c2 , a skyrmion lattice is predicted to have a 
lower energy than a vortex lattice Q • The currently most 
convincing case for a p-wave superconductor is Sr2Ru04 
[l3| . another candidate is UGe 2 fhH ]. 

Our central result is the prediction that one can dis- 
tinguish a skyrmion lattice from a vortex lattice by con- 
sidering the melting curve of the lattice. Our results for 
the latter are summarized in the phase diagram for clean 
skyrmion flux lattices shown in Fig. [lja). This is strik- 



ingly different from that for the vortex case, Fig. QJb): 
The skyrmion lattice never melts near H c \, while vortex 
lattices always melt near H c i. 
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FIG. 1: (a) External field (H) vs. temperature (T) phase 
diagram for skyrmion flux lattices. There is a direct transition 
from the skyrmion flux lattice to the Meissner phase. The 
theory predicts the shape of the melting curve only close to 
jT c ; the rest of the curve is an educated guess, (b) Same as (a) 
for vortex flux lattices. The vortex flux lattice always melts 
before H reaches H c i. 

The physics behind this result is simple. The inter- 
action between vortices falls off exponentially at large 
distances [HI]. As a result, when H — > H c i from above, 
which causes the vortex lattice constant R to grow, the 
elastic constants vanish exponentially with R as R — ► oo, 
which, in turn, causes the root-mean-square displacement 
fluctuation -\/|u(cc)| 2 to diverge exponentially. Conse- 
quently, -\/|ix(:e)| 2 ^>R&sR-^ooH~ > H c i. The Lin- 
demann criterion [lij then implies that the vortex lattice 
must melt with decreasing H before one reaches Hci ■ 

For skyrmions, as first shown numerically in Ref.U, and 
confirmed analytically below, the interaction potential 
falls off only as l/R. As a result, the Lame coefficients 
/i and Al of the skyrmion lattice vanish only as l/R? as 
R -> oo. This leads to ^/\u{x)\ 2 oc i? 3 / 4 < R. Hence, 
the skyrmion lattice never melts as H c i is approached. 
The shape of the phase diagram alone thus distinguishes 
skyrmion lattices from vortex lattices. 
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FIG. 2: Numerical data for the energy per skyrmion per unit 
length (circles) together with the best fit to a pure 1/R be- 
havior (dashed line) from Ref. |fj|, and the analytic solution 
given by Eq. {TJ (solid line). 



All of these results follow from an analytic solution of 
the Euler-Lagrange equations that give the minimum en- 
ergy configuration for a skyrmion lattice in the /3-phase of 
a p-wave superconductor. We have obtained an asymp- 
totic solution in the limit of large lattice spacing R ^> A, 
where A is the London penetration depth. Near H c \ the 
superconductor is always in this limit. As in the numeri- 
cal work of Ref. 0, we replace the hexagonal unit cell one 
expects in the skyrmion lattice by a circle of the same 
area. We expect this approximation to preserve the cor- 
rect scaling of the energy. 

For the first three terms of the asymptotic expansion 
of the energy E per skyrmion per unit length we obtain 

E/Eo = 2+(8%/6/3) (A/i?)-(16/3) (X/R) 2 ln(i?/A). (1) 

Here E = ($ /47rA) 2 , with $ = Trhc/e the flux quan- 
tum. This analytic result is in excellent agreement with 
the numerical solution reported in Ref. [4|, see Fig. [2] 
Notice that there are no free parameters in Eq. (JTJ) . 

In the remainder of this Letter we sketch the derivation 
of the above results. A complete account of the calcula- 
tions will be given elsewhere [13] ■ We begin by reviewing 
the formulation of the skyrmion lattice problem 0] . 

The spin part of the order parameter for a p-wave 
superconductor is a complex 3-vector ip{x). In a large 
regime of Landau theory parameter space (the so-called 
/3-phase), all low-energy configurations of ip(x) can be 
written in the form i/>(a;) = ipo(n(x) + irh(x)), where 
n(x) and rh(x) are real, mutually orthogonal unit vec- 
tors that vary slowly in space, and ?po is a constant 

In the ground state, n and rh are constants. Slow 



(diiy + (n • d t rh - a) 2 + (V x a) 2 
-2h ■ (V x a)] . (2) 



Here we use dimensionless units where distance, vector 
potential a, magnetic field h, and energy are measured in 
units of A, <E> /27rA, <& /27rA 2 , and $o/327r 3 A, respectively 
I = n x m, and the last term in the London energy, 
Eq. |2]), represents the coupling of the external magnetic 
field h to the magnetic induction b = V x a. 

A vortex is a low-energy configuration of h(x) and 
rh(x) in which I is a constant (i.e., in which h(x) and 
rh(x) span the same plane for all x), but h(x) and rh(x) 
rotate by 27rn (n integer) as one follows their evolution 
around any closed spatial path that encircles the path of 
the vortex core. Such a configuration necessarily has a 
singularity on the path of the core. 

The coupling between n(x), rh(x), and a makes the 
formation of vortices energetically favorable when h = 
(0, 0, h) is sufficiently strong. In conventional type-II su- 
perconductors, this leads to the appearance of a spon- 
taneous vortex flux lattice, where b is parallel to h, and 
n{x), rh(x), and b{x) all become spatially periodic func- 
tions of x and y (x = (x,y,z)). This two-dimensional 
spatially periodic structure forms a hexagonal lattice. 

For s-wave superconductors, vortices are the only pos- 
sible topological defects. For p-wave superconductors, 
the 3-d freedom of the vectors h{x) and rh(x) also al- 
lows for skyrmions, which are non-singular configurations 
of n(x) and rh(x). They look like n = 2 vortices at the 
unit cell boundary. Inside the cell, h(x) and rh{x) move 
out of their common plane, so l(x) is no longer constant. 
In the simplest case the skyrmion is cylindrically sym- 
metric, and l{x) forms an angle 9(x) with some central 
axis (which we will take to be the z-axis). While 0(x) 
changes from to 7r as one moves from infinity to the 
skyrmion axis, n(x) and rh(x) rotate around I by 47r on 
any loop enclosing the skyrmion axis. 

This picture leads to the characterization of a cylindri- 
cally symmetric skyrmion in terms of the single, as yet 
undetermined, angle 9(x). In polar coordinates (r, (p) the 
vector fields ii(x), m(x) and l(x) are given by 

I = e z cos 9{r) + e r sin 8(r), 
h = (e z sin 6(r) — e r cos 6(r)) sin <p + e v cost/?, 
rh = (e z sin#(r) — e r cos0(r)) cos^ — e v simp. (3) 

If the field 0(x) minimizes the energy of this configura- 
tion, then that configuration is a local minimum of the 
London energy, Eq. ((2]) [J]. 

Inserting Eq. into Eq. yields the energy per 
unit length, in units of Eq, of a cylindrically symmetric 



3 



skyrmion in a region of radius R, 



E = - 



dr \ 



dr i 



9'(r)) 2 + ^ sin 2 0(r) 
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cos 9(r)) + a(r) 



dr r [a(r)/r + a'(r)]' 



(4) 



The three terms represent the energy of a nonmagnetic 
skyrmion, the supercurrent energy, and the magnetic en- 
ergy, respectively. Skyrmions in a lattice are not cylindri- 
cally symmetric, since the lattice is not. However, since 
a hexagon is well approximated by a circle of the same 
area, we follow Ref. 4j by so approximating the unit cell. 

The Euler-Lagrange equations that result from mini- 
mizing the above energy read 
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a"(r) + -a'{r)-^a(r) 



a(r) + - [l + cos0(r)l . (5b) 
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This set of coupled, nonlinear ODEs must be solved 
subject to the boundary condition 9(r = 0) = 7r, and 
0(r = R) = 0. In Ref. |4j this was done numerically using 
finite elements methods. Here we show that for large R 
(i.e., near the lower critical field -ff c i) an analytic solution 
can be obtained in the form of an asymptotic expansion 
in powers of 1/R. To zeroth order, i.e, for R — > oo, the 
l.h.s. of Eq. (|5b|) vanishes, so a is given by 

o«,(r) = - [1 + cos 9(r)]/r. (6a) 

Inserting this in Eq. (|5a|) , the resulting ODE is solved by 

0oo(r) = 2arctan(^/r), (6b) 

for any £, All of these solutions fulfill the boundary con- 
dition, and £ is undetermined at this point. 

We write 0(r) = 9 00 {f) + 86(r) and a(r) = aoo(r) + 
Sa(r), and require |£0(r)| -C 1 and |<5a(r , )| -C |1 + 
cos0oo(r)|/r. We do not require \58\ <C 0oo, which is cru- 
cial for the success of our perturbative method. Inserting 
this into the ODEs, one sees that 59 can be written 

59{r) = (l/f)g(r/£) + (l/£ 4 )h(r/£) + 0(l/£ 6 ), (7) 

while 5a can be expressed in terms of 0OO, 59, and its 
derivatives. It will turn out that £ oc vR, so this is the 
desired expansion in powers of 1/R. For the function g 
we obtain a linear, inhomogeneous ODE of second order, 



g "(u) + -g'(u)- U 
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(8) 



The physical solution is the one that diverges linearly for 
large arguments. Standard methods give 

, , -4 u[u 2 (4 + u 2 ) + 2(1 + u 2 )\n(l + u 2 )} 

s(«) = T jYTut) 2 • (9) 



The parameter £ is now determined by the require- 
ment 0(r = R) = 0, which yields, to this order, £ = 
(2/3) 1 / 4 R 1 / 2 . To the same order, 5a is just given by the 
l.h.s. of Eq. (|5b|) with a replaced by a^, 



(10) 



5a{r) = (16r/€ 4 )(l + r 2 /£ 2 )- 3 + 0(l/£ 5 ). 



Inserting these results in Eq. (j4| and performing the 
integrals yields the first two terms on the r.h.s. of Eq. 
©. This method can be continued order by order. At 
the next order it yields the final term in Eq. (fT]). 

From the energy per skyrmion, Eq. (JTJ) , we can calcu- 
late the external field dependence of the equilibrium lat- 
tice constant, R(H). This is done by minimizing the en- 
ergy per unit volume, which is the energy per unit length 
per skyrmion, Eq. {T]), divided by the area per skyrmion, 
ttR 2 , minus an energy density gain of 2^qH/ttR 2 due 
to the external field. The latter is obtained from the 
h ■ (V x a) term in Eq. ([2]) by noting that the magnetic 
flux J dxdy z ■ (V X a) = 2<&o for each skyrmion in the 
lattice. This yields the Gibbs free energy density 



g(R) = {K/4tt 2 ) \-5/R 2 +4V6A/(3i? 3 ) 



(11) 



where K = $§/2ttA 2 and <5 = H/H cl - 1 with H cl = 
K/2<S> . 

For H < H c i, g(R) is minimized by R — > oo; i.e., there 
are no skyrmions. For H > H c \, on the other hand, the 
energy density is minimized at a finite R given by 



R = Ro = 2V6X/5. 



(12) 



Thus, H c i is the lower critical field at which the skyrmion 
lattice first forms. This implies for the spatially averaged 
magnetic induction, which is the flux per skyrmion 2$o 
divided by the unit cell area, 



B(H) = 2$„/K 2 ) = 5 2 H cl /3. 



(13) 



Hence, B(H) is horizontal near H c \, while it is vertical 
for vortex lattices. This result, with a slightly different 
numerical prefactor, was first obtained numerically in Q]. 

We now calculate the elastic properties of the skyrmion 
lattice. By symmetry, the elastic Hamiltonian of a hexag- 
onal lattice of lines parallel to the z-axis in 3-d is 

Hq1 = \ J dX ( 2 ^ UQ < 9 " a < 9 + ^L M "a u /3/3 + K tilt \d z u\ 2 ) , 

(14) 

where u a p = (dpu a + d a up) /2 is the strain tensor, a, (3 £ 
{x,y}, and u only has x and y components. Al, and 
-Ktiit are the shear, bulk, and tilt moduli. 

Consider the energy change due to a dilation of the 
lattice, R — > R£ = R (l + e) with e <C 1. This corre- 

and a strain ten- 



sponds to a displacement u(x) = e x j 
lar to z. Inserting this in fl" e i gives E^u/V = 2(/U+Al) e 2 , 



M>Y //,, , = eJ^fl, with x± the projection of a; perpendicu- 
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and comparing with the dilation energy implied by Eq. 

(ED, 

^ = g (RZ) - g(R ) = \ g"(R Q ) (eR ) 2 = e 2 . 

gives 

/i + Al = K 5 3 / I92n 2 X 2 . (15) 

Obtaining fi is more difficult since we have approxi- 
mated the unit cell by a circle. We use a heuristic ap- 
proach. Equation (fTTj) is of the form that would result if 
the skyrmions interacted via a nearest-neighbor pair po- 
tential U p (r) oc KX/r. Pretending that U p {r) is the ori- 
gin of the skyrmion energy (which should give the correct 
scaling of /i with h) makes it straightforward to calculate 
H, and comparing the result with Eq. (fT5|) yields Al- We 
obtain 

H~ X L = KS 3 /X 2 xO(l). (16) 

-Ktilt w e obtain by considering a uniform tilt of the 
skyrmion axes away from the z-axis by a small angle 
6 = 1 9*1*1 . The tilt energy in Eq. (fl~4"l) is the change in 
the b ■ h energy in Eq. @. Since, per skyrmion per unit 
length, this energy is (in ordinary units) — &qH cos 0/2it, 
tilting changes this energy by &qH(1 — cos#)/27r w 
$ H6 2 /4n = $ H\d z u\ 2 /4tt. Dividing this result by the 
unit cell area, using Eq. (TTZ]) for Rq, and identifying the 
result with the if t ilt term in Eq. fT4"|) gives if t ilt near 

K m = {H cl 5) 2 /I2ir. (17) 

We can now calculate the mean-square positional fluc- 
tuations (|w(a;)| 2 ) by Fourier transforming Eq. ([14")) and 
using the equipartition theorem. This gives 

(|„(x)| 2 ) T = ^ V(M<zi+^W) (18) 

qGBZ 

for the transverse fluctuations, and the same expression 
with fx — > (2fi + Al) for the longitudinal ones. q± and q z 
are the projections of q orthogonal to and along z, respec- 
tively. The Brillouin zone (BZ) of the skyrmion lattice 
is a hexagon of edge length oc 1/Rq in the plane perpen- 
dicular to z, and infinitely extended in the z-direction. 

Since /i and Al are comparable in magnitude, the same 
is true for the longitudinal and the transverse contribu- 
tions to (|tt(a;)| 2 ). Performing the integral over q z yields 

(|«(x)| a > ~ (|«(«)| a ) T = / |% , k f ■ (19) 

7bz 8tt z Vjtm tilt £1 

A change of variables, q± = w/Rq, and using Eqs. (|12p. 
(fTB"]) . and (fl~7|) gives, as claimed in the Introduction, 

(l-it^)! 2 ) = 0(1) x ksT/y/iiKtnt Ro 

= 0(1) x kBT/XH^S 3 /' 2 oc R 3/2 . (20) 



The Lindemann criterion for melting is Tl = 
(|zx(a;)| 2 )/i? 2 > T c = 0(1). In our case, 

T L = k B T5 1/2 /H 2 cl X b/2 x 0(1). (21) 

We see that, as claimed in the Introduction, the skyrmion 
lattice never melts as H — > H c \ from above (i.e., as 6 — > 
0), since the Lindemann ratio vanishes in that limit. 

We finally determine the shape of the melting curve 
H m (T) near the superconducting T c . Since, in mean field 
theory, H cl oc (T c - T), and A oc l/y/T c -T 0, we find 
from Eq. (jUJ by putting T L = O(l), 

Hm — H c i oc (T c — T) 5//2 . (22) 

The melting curve thus quickly rises above H c \ with de- 
creasing temperature, as shown qualitatively in Fig.flja). 
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